*     function out of Litkouhi and Beck (1982)
*     RF..............L............B.......82(r1,r2)=H(r1*r2,1/r2)=H(X,p)
*     equations (20) and (15), (16), (17a)
      DOUBLE PRECISION FUNCTION RFLB82(R1,R2)
      IMPLICIT NONE
      DOUBLE PRECISION R1, R2

      DOUBLE PRECISION RPI
      PARAMETER (RPI=3.1415926535)

      DOUBLE PRECISION RX, R3, R4, RH1, RH2, RH3, RP, RP1, RP3,
     &                 RAQUAD, RBQUAD, RCQUAD
      DOUBLE PRECISION RFERFC, RFLBSU

      RX=ABS(R1*R2)
      R3=ABS(R2)
      IF (R3.LT.1.E-8) THEN
         RFLB82=0.
      ELSE
         RP=1./ABS(R2)
         IF (RX.GT.3.) THEN
*           (equation 20)
            RFLB82=RFERFC(RX*RP)
         ELSE IF (RX.LT.1. AND. (3.*ABS(RP-1.)-RX).LT..3) THEN
*           interpolate between (rx;rx/3+1) (rx;1) and (rx;1-rx/3)
*           H(rx,rx/3+1.1)
            RP1=MAX(RX/3.+1.1, 9.-RX/3.)
*           (equation 15)
            R3=RX*RX*RP1*RP1
            R4=1./RP1
            RH1=2./RPI*RFLBSU(R3,R4)
*           H(rx,1)   (equation 17a)
            RH2=.5-.5*(1.-RFERFC(RX))**2
*           H(rx,.9-rx/3)
            RP3=MIN(RX/3.+1.1, .9-RX/3.)
*           (equation 16)
            R3=RX*RX
            R4=RP3*RX
            RH3=1.-(1.-RFERFC(RX))*(1.-RFERFC(R4))-
     &         2./RPI*RFLBSU(R3,RP3)
            IF (ABS(RP-RP1).LT..0001) THEN
               RFLB82=RH1
            ELSE IF (ABS(RP-1.).LT..0001) THEN
               RFLB82=RH2
            ELSE IF (ABS(RP-RP3).LT..0001) THEN
               RFLB82=RH3
            ELSE
               RAQUAD=((RH1-RH2)/(RP1-1.)-(RH2-RH3)/(1.-RP3))/(RP1-RP3)
               RBQUAD=(RH1-RH2)/(RP1-1.)-RAQUAD*(RP1+1.)
               RCQUAD=RH1-(RBQUAD+RAQUAD*RP1)*RP1
               RFLB82=(RAQUAD*RP+RBQUAD)*RP+RCQUAD
            END IF
         ELSE
            IF (RP.GT.1.001) THEN
               IF (RP.GT.1000.) THEN
                  RFLB82=0.
               ELSE IF (RX.LT..0001) THEN
                  RFLB82=2./RPI*ATAN(1/RP)
               ELSE IF (RX.GT.1.) THEN
*                 own initiative (because of exploding values)
                  R3=RX*RP
                  R4=RFERFC(R3)
                  RFLB82=R4*(1.-.5*R4)
               ELSE IF (-RX+3*RP.LT.3.) THEN
*                 own initiative (because of exploding values)
                  R3=RX*RP
                  R4=RFERFC(R3)
                  RFLB82=R4*(1.-.5*R4)
               ELSE
*                 (equation 15)
                  R3=RX*RX*RP*RP
                  R4=1./RP
                  RFLB82=2./RPI*RFLBSU(R3,R4)
               END IF
            ELSE IF (RP.LT..999) THEN
               IF (RP.LT..0001) THEN
                  RFLB82=1.
               ELSE IF (RX.LT..0001) THEN
                  RFLB82=1.-2./RPI*ATAN(RP)
               ELSE IF (RX+3*RP.GT.3.) THEN
*                 own initiative (because of negative values)
                  R3=RX*RP
                  R4=RFERFC(R3)
                  RFLB82=R4
               ELSE
*                 (equation 16)
                  R3=RX*RX
                  R4=RP*RX
                  RFLB82=1.-(1.-RFERFC(RX))*(1.-RFERFC(R4))-
     &               2./RPI*RFLBSU(R3,RP)
               END IF
            ELSE
*              RP about 1 (equation 17a)
               RFLB82=.5-.5*(1.-RFERFC(RX))**2
            END IF
         END IF
      END IF
*     H(X,-p)=-H(X,p)  (equation 18b)
*     H(-X,p)=+H(X,p)  (equation 18a)
*     --> H(-r1*r2,1./r2)=H(r1*r2,1./r2)
*     --> H(r1*(-r2),1./(-r2))=H(r1*r2,1./(-r2))=-H(r1*r2,1.)
      IF (R2.LT.0.) RFLB82=-RFLB82
      END
